Spectral fluctuations of billiards with mixed dynamics: from time series to 

superstatistics 
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A statistical analysis of the eigenfrequencies of two sets of superconducting microwave billiards, 
one with mushroom-like shape and the other from the familiy of the Limacon billiards, is presented. 
These billiards have mixed regular-chaotic dynamics but different structures in their classical phase 
spaces. The spectrum of each billiard is represented as a time series where the level order plays the 
role of time. Two most important findings follow from the time-series analysis. First, the spectra 
can be characterized by two distinct relaxation lengths. This is a prerequisite for the validity of 
the superstatistical approach which is based on the folding of two distribution functions. Second, 
the shape of the resulting probability density function of the so-called superstatistical parameter 
is reasonably approximated by an inverse \ 2 distribution. This distribution is used to compute 
nearest-neighbor spacing distributions and compare them with those of the resonance frequencies 
of billiards with mixed dynamics within the framework of superstatistics. The obtained spacing 
distribution is found to present a good description of the experimental ones and is of the same 
or even better quality as a number of other spacing distributions, including the one from Berry 
and Robnik. However, in contrast to other approaches towards a theoretical description of spectral 
properties of systems with mixed dynamics, superstatistics also provides a description of properties 
of the eigenfunctions. Indeed, the inverse x 2 parameter distribution is found suitable for the analysis 
of experimental resonance strengths in the Limagon billiards within the framework of superstatistics. 

PACS numbers: 02.50.-r, 05.40.-a, 05.4S.Mt, 05.45.Tp, 03.65.-w 



I. INTRODUCTION 



Integrable Hamiltonian dynamics is characterized by 
the existence of as many conserved quantities as de- 
grees of freedom. Each trajectory evolves on an invari- 
ant hyper-torus in the phase space [HQ- In contrast, 
chaotic systems are ergodic; almost all orbits fill the en- 
ergy shell in a uniform way. Physical systems with inte- 
grable and fully chaotic dynamics, respectively, are, how- 
ever, exceptional. A typical Hamiltonian system shows 
a phase space in which regions of regular motion and 
chaotic dynamics coexist. These systems are known as 
mixed systems. Their dynamical behavior is by no means 
universal, as is the case for fully regular and fully chaotic 
systems. If we perturb an integrable system, most of 
the periodic orbits on tori with rational frequencies dis- 
appear. However, some of these orbits persist. Elliptic 
periodic orbits appear surrounded by islands. They cor- 
respond to librational motions around these periodic or- 
bits and reflect their stability. The Kolmogorov-Arnold- 
Moser (KAM) theorem states that invariant tori with a 
sufficiently incommensurate frequency vector are stable 
with respect to small perturbations. Numerical simula- 
tions show that when the perturbation increases more 
and more tori are destroyed. For large enough perturba- 
tions, there are locally no tori in the considered region 
of phase space. The break-up of invariant tori leads to a 
loss of stability of the system, that is, to chaos. There are 
three main scenaria of transition to global chaos in finite- 
dimensional (nonextended) dynamical systems, one via a 



cascade of period-doubling bifurcations, a Lorenz system- 
like transition via Hopf and Shil'nikov bifurcations, and 
the transition to chaos via intermittency 3, 4]. It is nat- 
ural to expect that there are other (presumably many 
more) such scenaria in extended (infinite-dimensional) 
dynamical systems. 

In quantum mechanics, the specification of a wave 
function is always related to a certain basis. In integrable 
systems the eigenbasis of the Hamiltonian is known in 
principle. In this basis, each eigenfunction has just one 
component. That obviously indicates the absence of com- 
plexity. In the nearly ordered regime, mixing of quantum 
states belonging to adjacent levels can be ignored and the 
energy levels are uncorrelated. The level-spacing distri- 
bution is well described by that for random numbers gen- 
erated by a Poissonian process, exp(— s), where s is the 
spacing between adjacent energy levels rescaled to unit 
mean spacing D. For a Hamiltonian with a chaotic clas- 
sical limit, on the other hand the wavefunction compo- 
nents are on average uniformly distributed over the whole 
basis. Berry [B| conjectured that the wavefunctions of 
chaotic quantum systems can be represented as a formal 
sum over elementary solutions of the Laplace equation 
in which the real and imaginary parts of the coefficients 
are independent identically-distributed Gaussian random 
variables with zero mean and variance computed from the 
normalization. Bohigas et al. [6| put forward a conjec- 
ture (strongly supported by accumulated numerical and 
experimental evidence and also byrecent advances to- 
wards a proof of this conjecture |7j]) that the spectral 
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statistics of chaotic systems follow random-matrix theory 
(RMT, see [H HI). The properties of a chaotic Hamilton 
operator can thus be modeled by an ensemble of random 
Hermitian matrices H that belongs to one of three uni- 
versality classes, either the orthogonal, the unitary or the 
symplectic one and is called Gaussian orthogonal (GOE), 
unitary (GUE) and symplectic (GSE) ensemble, respec- 
tively. The theory is based on two main assumptions: the 
matrix elements are independent identically-distributed 
random variables and their distribution is invariant un- 
der unitary transformations. This leads to the Gaussian 
probability density distribution for the matrix elements 

P(H) = ^-e X p[-r l Tr(H^H)], (1) 

where Z(rj) = J exp [— r/Ti (ii^_ff)] dr\ is the normaliza- 
tion constant. The Gaussian distribution is also obtained 
by maximizing the Shannon entropy with the constraints 
of normalization and existence of the expectation value 
of Ti-(WH), see, e.g. @, QjJ. Information about the 
statistical properties of the eigenvalues and/or eigenvec- 
tors of the matrix H can be obtained by integrating over 
the undesired variables. There is strong evidence by now, 
that indeed the spectral correlation functions of a chaotic 
system are well decribed by those obtained from Eq. (TT|) 
and determined solely by the global symmetries of the 
system such as time-reversal invariance and the value 
of the spin. Among the measures representing spectral 
correlations, the nearest-neighbor level-spacing distribu- 
tion (NNSD) p(s) has been studied extensively so far. 
For the random matrix ensembles Eq. |T]) it is well ap- 
proximated by the Wigner-Dyson distribution, namely 
Pp(s) — aps 13 ex.p(—bps 2 ), where (3 (= 1, 2, and 4 for 
the orthogonal, the unitary, and the symplectic ensem- 
bles, respectively) characterizes the universality classes. 
The coefficients ap and bp are determined by the nor- 
malization conditions J pp(s)ds — j sp@(s)ds = 1, as 
a x = tt/2, a 2 = 32/tt 2 , a 4 = tt 3 2 18 /3 6 , b x = vr/4, b 2 = 
4/7r,and 64 = 64/97T. For s«l, the distribution function 
is proportional to s@ , which implies that adjacent energy 
levels repel each other. This behavior may be attributed 
to the mixing between the two states related with these 
levels. 

So far in the literature, there is no rigorous statistical 
description for the transition from integrability to chaos. 
The nature of the stochastic transition is more obscure 
in quantum than in classical mechanics, as the assump- 
tions that lead to the RMT description do not apply to 
mixed systems. The Hamiltonian of a typical mixed sys- 
tem can be described as a random matrix where some of 
its elements are randomly distributed and some of them 
might be non-random. Moreover, the matrix elements 
need not all have the same distributions and may or may 
not be correlated. Thus, the RMT approach is a difficult 
route to follow. Comprehensive semiclassical computa- 
tions have been carried out for Hamiltonian quantum 
systems, which on the classical level have a mixed phase 
space dynamics (see, e.g. [ll[ and references therein). 



There have been several proposals for phenomenologi- 
cal random matrix theories that interpolate between the 
Wigner-Dyson RMT and banded random matrices with 
an almost Poissonian spectral statistics. The standard 
route for the derivation is to sacrifice basis invariance 
but keep matrix-element independence. The first work 
in this direction is due to Rosenzweig and Porter [I3 |. 
They model the Hamiltonian of a mixed system by a su- 
perposition of a diagonal matrix with random elements 
and a matrix drawn from a GOE. Accordingly, the vari- 
ances of the diagonal elements of the total Hamiltonian 
are not twice that of the off-diagonal ones, as in the GOE 
case. Hussein and Pato [l3T ] used the maximum entropy 
principle to construct such ensembles by imposing ad- 
ditional constraints. Also, ensembles of banded random 
matrices whose entries are equal to zero outside a band 
of width b along the principal diag onal have been used 
to model mixed systems [3 [IE EE El El ■ 

Another route for generalizing RMT is to conserve 
base invariance but allow for the correlation of matrix 
elements. This has been achieved by maximizing non- 
extensive entropies subject to the constraint of a fixed 
expectation value of Tr(H^H), see 0, H H [H [H, H 
1 2 51 . Recently, an equivalent approach was presented in 
[26l [27J , which is based on the method of superstatistics 
(statistics of a statistic) proposed by Beck and Cohen 
[28|. This formalism has been elaborated and applied 
successfull y to a wide variety of physical problems, e.g., 
in [H, [MiM H [H, HI HE HI. In thermostatics, su- 
perstatistics arises from weighted averages of ordinary 
statistics (the Boltzmann factor) due to fluctuations of 
one or more intensive parameters (e.g. the inverse tem- 
perature). Its application to RMT assumes the spectrum 
of a mixed system as made up of many smaller cells that 
are temporarily in a chaotic phase. Each cell is large 
enough to obey the statistical requirements of RMT but 
is associated with a different distribution of the parame- 
ter 77 n Eq. |T]) according to a probability density f(rj). 
Consequently, the superstatistical random-matrix ensem- 
ble used for the description of a mixed system consists of 
a superposition of Gaussian ensembles. Its joint probabil- 
ity density distribution of the matrix elements is obtained 
by integrating the distribution given in Eq. |T]) over all 
positive values of 7/ with a statistical weight f(rf), 

f°° , , exp [-rjTr (WH)} 

p{h) = J f( v ) yi ' 'Mr,- (2) 

Despite the fact that it is hard to make this picture rig- 
orous, there is indeed a representation which comes close 
to this idea [HEI. 

The present paper is concerned with a justification for 
the use of the above-mentioned superstatistical general- 
ization of RMT in the study of mixed systems, based 
on the representation of their energy spectra in the form 
of discrete time series in which the level order plays the 
role of time. The representation of the suitably trans- 
formed eigenvalues of a quantum system as a time series 
has recently allowed to determine the degree of chaotic- 
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ity of the dynamics of the system [42l. l43l lij. Hrj]. 
We have thus been motivated by the work of Beck, Co- 
hen and Swinney [471 ] concerning the derivation of su- 
perstatistics starting from time-series. Superstatistical 
thermostatics results as a convolution of two statistics, 
one characterized by the Boltzmann factor and the other 
corresponding to inverse-temperature fluctuations. This 
requires the existence of two relaxation times. We ap- 
ply the arguments of [47| by representing the spectra of 
mixed systems as discrete time series in which the role 
of time is played by the level ordering. In Section II, 
we consider two billiards with mushroom-shaped bound- 
aries as representatives of systems with mixed regular- 
chaotic dynamics and three with the shape of Limacon 
billiards, one of them of chaotic and two of mixed dynam- 
ics. The quantum eigenvalues and statistical properties 
of the eigenfunctions were obtained experimentally by 
exploiting the equivalence of the Schrodinger equation 
of a plane quantum billiard and the Helmholtz equa- 
tion for the electric field strength in a cylindrical mi- 
crowave resonator for wave lengths longer than twice the 
height of the resonator. The billiards with mixed dynam- 
ics have classical phase spaces of different structures for 
the two families of billards. The "time-series" analysis 
of their spectra manifests the existence of two relaxation 
lengths, a short one defined as the average length over 
which energy fluctuations are correlated, and a long one 
that characterizes the typical linear size of the hetero- 
geneous domains of the total spectrum. It is performed 
in an attempt to clarify the physical origin of the het- 
erogeneity of the matrix-element space, which justifies 
the superstatistical approach to RMT. The second main 
result of this section is to derive a parameter distribu- 
tion f{rf), which is introduced in Eq. (0). This paves 
the way for the generalization of the Wigner surmise to 
superstatistics concerning the nearest-neighbor spacing 
distribution (NNSD). We then apply the deduced gener- 
alized Wigner surmise in a phenomenological analysis of 
the NNSD to the measured resonance frequencies of the 
microwave resonators. Section IV introduces supersta- 
tistical generalizations for the Porter-Thomas distribu- 
tion of partial widths. Then, the corresponding formulas 
for the resonance strengths are used to analyze the ex- 
perimental resonance-strength distributions in the mixed 
Limagon billiards. A brief summary of the main results 
is given in section V. 



II. MUSHROOM AND LIMAgON BILLIARDS 

Billiards can be used as simple models in the study 
of Hamiltonian systems. They consist of a point parti- 
cle which is confined to a container of some shape and 
reflected elastically on impact with the boundary. The 
shape determines whether the dynamics inside the bil- 
liard is regular, chaotic or mixed. The best-known ex- 
amples of chaotic billiards are the Sinai billiard (a square 
table with a circular barrier at its center) and the Buni- 



movich stadium (a rectangle with two circular caps) [48j |. 
Neighboring parallel orbits diverge when they collide 
with dispersing components of the billiard boundary. In 
chaotic focusing billiards, neighboring parallel orbits con- 
verge at first, but divergence prevails over convergence on 
average. Divergence and convergence are balanced in in- 
tegrable billiards such as circles and ellipses. 

Recently Bunimovich introduced the so-called 'mush- 
room' billiard [49] with the novel feature of a well- 
understood divided phase-space comprising a single in- 
tegrable region and a single ergodic one. We restrict 
ourselves here to mushroom billiards which consist of a 
semicircular region, the 'hat' and a 'stem', which is sym- 
metrically attached to its base. As the width of the stem 
varies from zero to the diameter of the hat, there is a 
continuous transition from integrability (the semicircle 
billiard) to ergodicity (in case of a rectangul stem the sta- 
dium billiard). In mushroom billiards, the regular region 
has a well-defined semicircular border about the center of 
the hat with radius equal to half the width of the stem. 
It is composed of those trajectories in the hat that never 
cross this border and therefore remain in the hat forever. 
Their integrability is due to the conservation of the reflec- 
tion angle for collisions with the semicircular boundary. 
The chaotic component consists of trajectories that en- 
ter the stem of the mushroom billiard. In contrast to 
most other mixed systems, the dynamics of mushroom 
billiards is free of the usual hierarchies of KAM islands 
about integrable islands in phase space. Because of its 
sharply-divided phase space, mushroom billiards can be 
thought of as an ideal model for the understanding of 
mixed dynamics. They indeed have already been under 
active research (5(J HH, [H, HI] . 

The Limagon billiard is a closed billiard whose bound- 
ary is defined by the quadratic conformal map of the unit 
circle z to w, 

w = z + Xz 2 , \z\ = 1. (3) 

The shape of the billiard is controlled by a single pa- 
rameter A with A = corresponding to the circle and 
A = 1/2 to the cardioid billiard [54} . For < A < 1/4. the 
Limagon billiard has a continuous and convex boundary 
with a strictly positive curvature and a collection of caus- 
tics near the boundary (5f| [5(|. At A = 1/4, the bound- 
ary has zero curvature at its point of intersection with 
the negative real axis, which turns into a discontinuity for 
A > 1/4. Accordingly, there the caustics no longer per- 
sist [11]. The classical dynamics of this system and the 
corresponding quantum billiard have been extensively in- 
vestigated by Robnik and collaborators [H, H3- They 
concluded that the dynamics in the Limagon billiard un- 
dergoes a smooth transition from integrable motion at 
A = via a soft chaos KAM regime for < A < 1/4 to a 
strongly chaotic dynamics for A = 1/2. 

Both families of systems have been studied experimen- 
tally in the quantum limit exploiting the analogy between 
a quantum billiard and a flat cylindric microwave billiard 
[53I . [58 . 59] . The electromagnetic resonances in a flat mi- 
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crowave cavity can directly be associated with quantum 
states in a quantum billiard of the same geometry. For 
the evaluation of statistical measures, a sufficiently large 
number of resonances is needed. Experimentally, this is 
not trivial since each resonance has a finite width, and 
as the level density increases with frequency, single res- 
onances can only be resolved up to a certain frequency. 
Hence, to measure as many states as possible, one has to 
reduce the width of the resonances by reducing the loss 
mechanisms. This is achieved by the use of superconduct- 
ing cavities, which are cooled down to the temperature 
of liquid helium, T = 4.2 K, in a bath cryostat. The 
cavities are either made of niobium or of copper plates, 
which have been galvanically covered with a layer of lead, 
whose thickness is several penetration depths of the elec- 
tromagnetic field. 

Electric field oscillations in the interior of the cavity 
can be excited via antennae. Using a vectorial network 
analyzer, the complex amplitude ratio of input to out- 
put signal from the cavity can be measured. Peaks in 
the modulus of the amplitude are found at resonance 
frequencies, corresponding to eigenmodes of the system. 
There, the wave field forms a standing wave, ft can be 
only excited when the antennae used in the measurement 
process are not near nodal lines of the wave function, 
where the amplitude is very small. Thus, usually several 
antennae are used to measure different spectra, enabling 
to identify all resonances up to a certain frequency. For 
flat resonators of a length scale of 30 cm, one collects 
approximately the first 800 eigenvalues with a very high 
precision. Experimental data of both kinds of systems 
considered in this paper have been obtained via this pro- 
cedure [6(3]. 

Two mushroom billiards have been recently investi- 
gated experimentally. In order to avoid symmetry ef- 
fects and bouncing ball orbits between parallel walls in 
the stem, their shape is of a half mushroom with a slant 
stem (see inset of Fig. [T] The ratio of the width of the 
stem to the diameter of the hat is 1:3 (2:3) for the small 
(large) mushroom billiard. The degree of chaos, which 
is the measure of all chaotic parts of the phase space, 
is 45.5 % (82.9 %), and the first 780 (938) resonances 
could be detected. In Fig. [T] we show a part of the 
spectra of the small and large mushroom billiard, respec- 
tively. Details of the experiment with the larger mush- 
room billiard can be found in [53j . Three desymmetrized 
cavities with the shape of billiards from the family of 
Limacon billiards have been constructed for the values 
A = 0.125,0.150,0.300 and the first 1163, 1173 and 942 
eigenvalues were measured, respectively. More details 
on these experiments are given in [58L 153 ] . To compare 
the statistical properties of the eigenvalues with univer- 
sal predictions considered in the present paper, they have 
to be rescaled to unit mean spacing. This is done by an 
unfolding procedure using Weyl's formula [6l|, which re- 
lates the billiard area and circumference to the number 
of resonance frequencies below a given one. 

As outlined above, the phase space structure of the 
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FIG. 1: Part of the transmission spectrum of the small (left 
panel) and large (right panel) mushroom billiard. The insets 
show the geometry of the billiards. 

billiards with mixed dynamics are different for the two 
families under consideration. The experimental NNSD of 
the larger of the two mushroom billiards exhibits [53j a 
statistically significant dip at s w 0.7. It vanishes when 
the contribution of the two shortest periodic orbits is 
subtracted. Such a dip has never been observed in the 
spectra of other billiards with mixed dynamics, including 
the Limagon billiards considered here. We shall show 
below that nevertheless the statistical properties of the 
Limagon billiards are indistinguishable from those of the 
mushroom billiards after removal of the contribution of 
these periodic orbits, in spite of the difference of their 
phase space structure. 

III. TIME-SERIES REPRESENTATION 

In this section the time series method used for the 
study of the fluctuations of the resonance spectra of the 
mushroom billiards is introduced. Representing energy 
levels of a quantum system as a discrete time series has 
been probed in a number of recent publications. Relano 
et al. [42] considered a sequence of energy levels as a 
discrete time series in which the energy played the role 
of time. They conjectured that the power spectra of 
chaotic quantum systems are characterized by 1/ / noise, 
whereas integrable quantum systems exhibit l// 2 noise. 
This conjecture was supported by numerical experiments 
which involved classical random-matrix ensembles and 
atomic nuclei. Moreover, the power spectrum of an ex- 
perimentally modeled quantum Sinai billiard exhibits a 
clear 1// noise through almost the whole frequency do- 
main. Mixed systems on the other hand, like the Limagon 
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billiard, the quartic coupled oscillator, and the kicked 
top, are characterized by a l/f a noise [H, |44|- In all 
these cases, the exponent a was related to the degree 
of chaos. Manimaran et al. [45[ recently developed a 
wavelet based approach to discrete time series and em- 
ployed it to characterize the scaling behavior of spectral 
fluctuations of random matrix ensembles, as well as com- 
plex atomic systems. Santhanam et al. [46| studied the 
spectra of atoms and Gaussian ensembles using the de- 
trended fluctuation analysis, which is a popular tool to 
study long range correlations in time series [62| . They 
showed that this analysis is related to the A3 statistics 
of RMT. 

Here we apply the time-series analysis to the study of 
the energy spectra of mixed systems from another point 
of view. Our goal is to test the hypothesis that quan- 
tum systems can be modeled by a generalization of RMT 
[HI, H3] based on the concept of superstatistics [28[ . We 
follow an approach recently proposed by Beck, Cohen 
and Swinney [47| . which describes how to proceed from 
a given experimental time series to a superstatistical de- 
scription. This approach allows one to check whether a 
time series contains two separate time scales, and also to 
extract the relevant probability densities of superstatis- 
tical parameters from the time series. 



A. Spectral relaxation lengths 

Superstatistical RMT assumes that the space of ma- 
trix elements consists of many spatial cells with differ- 
ent values of some intensive parameter, e.g. the inverse 
variance r\. In systems with mixed regular-chaotic dy- 
namics, the origin of this spectral heterogeneity is the 
possible partial conservation of an unknown or ignored 
symmetry. In the space of matrix elements each hetero- 
geneous domain comprises those matrix elements that 
couple states, which have similar properties with respect 
to this symmetry, where the typical size of the hetero- 
geneous cells T measures the correlation length in that 
space. The heterogeneity of the space of matrix elements 
presumably causes one in the structure of the spectrum. 
Each cell is assumed to reach local equilibrium very fast, 
i.e., the associated relaxation length r, which is defined 
as that length-scale over which energy fluctuations are 
correlated, is short. It may also be regarded as an op- 
erational definition for the average energy separation be- 
tween levels due to level repulsion. In the long-term run, 
the stationary distributions of this inhomogeneous sys- 
tem arise as a superposition of the "Boltzmann factors" 
of the standard RMT, i.e. e~ nTlH . The parameter 77 
is approximately constant in each cell for an eigenvalue 
interval of length T. In superstatistics this superposition 
is performed by weighting the stationary distribution of 
each cell with the probability density f(r)) to observe 
some value r\ in a randomly chosen cell and integrating 
over rj. Of course, a necessary condition for a supersta- 
tistical description to make sense is the condition t«T, 



because otherwise the system is not able to reach local 
equilibrium before the next change takes place. 

Our goal is to show that the behavior of a fictitious 
time series formed by the 780 and 938 resonances, re- 
spectively, in the two mushroom billiards, and the ap- 
proximately 1100 resonances in each of the three Limacon 
billiards is consistent with superstatistics. For this pur- 
pose the distribution f(rj) is derived by proceeding as in 
[471 ] . We extract the relaxation lengths (times) to local 
equilibrium r and the large length scale T on which the 
intensive parameter fluctuates and show that there is a 
clear scale separation of the spectral correlations in each 
billiard. 

First, let us determine the long time scale T. For this 
we divide the spacings series into N equal level-number 
intervals of size n. The total length of the spectrum is 
N ■ n. We then define the mean local kurtosis «(n) of a 
spacing interval of length n by 

N i=i (( s - s ) )i,N 

Here (• • • } i N = Y?k=(i-i)-n+i ' ' ' denotes a summation 
over an interval of length n starting at level spacing in, 
and s is either the local average spacing in each spac- 
ing interval or the global average s = 1 over the entire 
spacings series. We chose the latter one. In probability 
theory and statistics, kurtosis is a measure for the "flat- 
ness" of the probability distribution of a real- valued ran- 
dom variable. Higher kurtosis means that a larger part of 
the contributions to the variance is due to infrequent ex- 
treme deviations, as opposed to frequent modestly sized 
ones. A superposition of local Gaussians with local flat- 
ness three results in a kurtosis of three. We define the 
superstatistical level-number scale T by the condition 

k(T) = 3, (5) 

that is, we look for the simplest superstatistics, a su- 
perposition of local Gaussians [47| ■ If n is chosen such, 
that only one value of s is contained in each interval, 
then of course k(1) =1. If on the other hand n com- 
prises the entire spacing series, then we obtain the flat- 
ness of the distribution of the entire signal, which will 
be larger than 3, since superstatistical distributions are 
fat-tailed. Therefore, there exists a level-number scale T 
which solves Eq. ^ . Figure shows the dependence of 
the local flatness of a spacing interval on its length for the 
two mushroom and the three Limagon billiards. In the 
case of the chaotic Limagon billiard, in which A = 0.300, 
the quantity k does not cross the line of n = 3 for the 
considered values of n. It is expected that T — N in 
this case since the fluctuations in a chaotic (unfolded) 
spectrum are uniform. The values of T for the mixed 
billiards are given in Table IIII Al The short time scale, 
that is the relaxation time associated with each of the N 
intervals, was estimated in [471 ] from the small- argument 
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FIG. 2: The local kurtosis of the spacing intervals of the 
mushroom (left panel) and Limagon (right panel) billiards. 



TABLE I: Correlation lengths estimated from the time se- 
ries representing the frequency spectra of the mushroom and 
Limagon billiards. 
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exponential decay of the autocorrelation function 



of the time series s(t) under consideration. Figure [3] 
shows the behavior of the autocorrelation functions for 
the series of resonance-spacings of the two families of 
billiards. Quite frequently, the autocorrelation function 
shows single-exponential decays, C(n) = e~ n / T , where 
t > defines a relaxation " time" . A typical example is 
the velocity correlation of Brownian motion [63| . The au- 
tocorrelation functions studied here clearly do not follow 
this trend. For the systems with mixed dynamics, they 
decay rapidly from a value of C(0) = 1, change sign at 



FIG. 3: The autocorrelation of the spacings within spectral 
intervals for the mushroom (left panel) and Limagon (right 
panel) billiards. 

some n ~ 1 becoming negative, then asymptotically tend 
to zero. In an attempt to quantify the dependence of C s 
on n, we parametrized its empirical value in the form of 
a superposition of two exponentially decaying functions 

C s (n) = A ie - n/Tl + A 2 e- n/T2 (7) 

and (arbitrarily) fixed the superposition coefficient as 
A\ = 1.5 and A 2 = —0.5. The curves in Fig. [3] show 
the resulting parametrization. The best fit parameters 
are given in Table IIII Al We may estimate r as the 
mean values of t± and t 2 and conclude that r has a 
value slightly larger than 1 for each billiard. This is suf- 
ficient to conclude that the ratio T/r is large enough in 
each billiard to claim two well separated "time" scales in 
the level-spacings series, which justifies describing them 
within the framework of superstatistics. 

Both Figs. [2] and [3] as well as Table suggest that the 
evolution of the quantities n(n) and C s (n) along the time 
series for the two billiard families is the same despite the 
different behavior of their classical dynamics. 

B. Estimation of the parameter distribution 

The distribution f(rj) is determined by the spatiotem- 
poral dynamics of the entire system under consideration. 
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Beck et al. [47j have argued that typical experimental 
data are described by one of three superstatistical uni- 
versality classes, namely, \ 2 1 inverse \ 2 ■> or log-normal 
superstatistics. The first is the appropriate one if rj has 
contributions from v Gaussian random variables X\, . . 
. , X u due to various relevant degrees of freedom in the 
system. As mentioned before rj needs to be positive; this 
is achieved by squaring these Gaussian random variables. 
Hence, i] = J2i=i X? is x 2 distributed with degree v, 



L.9 



m = 



1 



I>/2) \2 Vo J 



\ v/2 



(8) 



The average of r\ is r) = J °° r\j(r\)dr\. The same consid- 
erations are applicable if ry -1 , rather than i], is the sum of 
several squared Gaussian random variables. The result- 
ing distribution f(rj) is the inverse x 2 distribution given 
by 



T{v/2) \ 2 J 



v/2 



n 



-I//2-2 p -vr, /2r, 



(9) 



where again rj is the average of 77. Instead of being a 
sum of many contributions, the random variable 77 may 
be generated by multiplicative random processes. Then 
In 77 = X)i=i m ^ is a sum of Gaussian random variables. 
Thus it is log-normally distributed, 



m = 



1 



[\n( V /v)] 2 /2v 2 



2ttv7] 



(10) 



which has an average fiy/w and variance /j, 2 w(w — 1), 
where w = exp(i; 2 ). 

Next, we need to determine which of these distribu- 
tions fits best that of the slowly varying stochastic pro- 
cess rj(t) described by the experimental data. Since the 
variance of superimposed local Gaussians (see remark af- 
ter Eq. ([5])) is given by ?7 _1 , we may determine the process 
r](t) from the series 



r?(i) 



1 



(11) 



Accordingly, the probability density f(rj) is determined 
from the histogram of the rj (i) values for all i\ the re- 
sulting experimental distributions are shown in Fig. [4] 
We compared them with the log-normal, the % 2 and the 
inverse x 2 distributions with the same mean (77) and vari- 
ance (?/ 2 ) — (v) 2 - The inverse \ 2 distribution fits the data 
significantly better than the other two distributions. 



IV. NEAREST-NEIGHBOR SPACING 
DISTRIBUTION 

This section focusses on the question whether the in- 
verse % 2 distribution of the superstatistical parameter 77 
in Eq. © is suitable for describing the NNSD of sys- 
tems in the transition out of chaos within the super- 
statistical approach to RMT. As mentioned above, the 



Mu^fcTonnn ftilliurd^- 



f.n-jJI bil!iard 




FIG. 4: Estimation of the parameter distribution for the su- 
perstatistical description of spectra of the mushroom (left 
panel) and Limagon (right panel) billiards. The solid lines 
represent the x 2 , the daotted lines the inverse \ 2 an d the 
dashed lines the log-normal distribution. 



NNSD of a chaotic system is well described by that of 
random matrices from the GOE, the Wigner surmise, if 
the system is chaotic and by that for Poisson statistics 
if it is integrable. Numerous interpolation formulas de- 
scribing the intermediate situation between integrability 
and chaos have been proposed [9( . One of the most pop- 
ular ones is that introduced by Brody [iH although it 
is purely phenomenological. This distribution coincides 
with the Wigner distribution for a fully chaotic and with 
Poisson's for an integrable system. It is known to pro- 
vide an excellent description for the NNSDs of numerous 
mixed systems. Another phenomenological distribution 
was proposed in [14j | and its usefulness was demonstrated 
for band random matrices. Lenz and Haake [65| derived a 
distribution based on the model of additive random ma- 
trices. Finally, Berry and Robnik elaborated a NNSD for 
mixed systems based on the assumption that semiclassi- 
cally the eigenfunctions and associated Wigner distribu- 
tions are localized either in classically regular or chaotic 
regions in phase space 66]. Accordingly, the sequences 
of eigenvalues connected with these regions are assumed 
to be statistically independent, and their mean spacing is 
determined by the invariant measure of the correspond- 
ing regions in phase space. The largest discrepancy be- 
tween the Brody and the Berry-Robnik distribution is 
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observed for level spacings s close to zero. While the for- 
mer vanishes for s = 0, the latter approaches a constant 
and nonvanishing value for s — > 0. In (53l . [59^ the exper- 
imental NNSDs for the measured resonance frequencies 
of the mushroom and the Limagon billiards, respectively, 
were compared to the Berry-Robnik (BR) distribution. 

It follows from Eq. ^ that the statistical measures of 
the eigenvalues of the superstatistical ensemble are ob- 
tained as an average of the corresponding ^-dependent 
ones of standard RMT weighted with the parameter dis- 
tribution f(p)- In particular, the superstatistical NNSD 
is given by [26| as 



p{s) 



fiv) Pw(v, s ) drj, 



(12) 



where pw{Vt s ) is the Wigner surmise for the Gaussian 
orthogonal ensemble with the mean spacing depending 
on the parameter 77, 



(13) 



For a x 2 distribution of the superstatistical parameter 77, 
one substitutes Eq. © into Eq. (|13[) and integrates over 
77. The resulting NNSD is given by 



p x *{v,s) 



(1+77 S 2 /^) 1+I V 2 ' 



(14) 



The parameter rj is fixed by requiring that the mean- 
level spacing (s) equals unity, yielding 



T 



V- 1 



r 



(15) 



For an inverse x 2 distribution of 77, given by Eq. ([9]) , one 
obtains the following superstatistical NNSD, 



Pi 



11 v X" 



2 Vo s 
T(u/2) 



r] us/2)" 2 K u/2 (y/ffou s ) 



(16) 

where K m (x) is a modified Bessel function [67J and r] 
again is determined by the requirement that the mean- 
level spacing (s) equals unity, 



47T 

Vo = — 



"3 



(17) 



Finally, if the parameter 77 has a normal distribution (10), 
then the NNSD 



Pio g n(w,s) 



2irv 



exp 



77s 



In' 



2tj 2 



(18) 

can only be evaluated numerically. 

We compared the resulting NNSDs given in Eqs. (|X4|>. 
((THJ) and [)18p with the experimental ones for the mush- 
room billiards and the two Limacon billiards with mixed 



1- -. 




11.6 



ti.a 



FIG. 5: Experimental NNSDs for the mushroom (left panel) 
and Limagon (right panel) billiards compared with the super- 
statistical distributions. The solid lines represent the x 2 j the 
dashed lines the inverse \ 2 y the dashed lines the log-normal 
and the short-dashed lines the Berry-Robnik distribution. 



dynamics. In Fig. [5] the experimental results are shown 
together with the superstatistical and the BR distribu- 
tions. The best fit values of the parameters ar e g iven in 
Table |TTJ We in addition used the entropy [68|, [63] 



S(P\p) 



P(x) In ^rldx 
p{x) 



(19) 



to express the difference between the theoretical Ansatz 
p{x) and the experimental distribution P(x). Such a ref- 
erence term has been discussed in the literature in the 
context of going from a discrete to a continuous system 
and is proportional to the limiting density of discrete 
points 70]. As is well known, the entropy S(P \p) is neg- 
ative everywhere except at its maximum, where it equals 
zero and P(x) = p(x). One can use this fact to find the 
value of the parameter of a theoretical distribution that 
has the least distance to the experimental distribution. 
The resulting best fit values of the parameters are given 
in Table ED Figure M and Tables HI] and [m] suggest the 
validity of the superstatistical distribution, especially for 
the nearly chaotic billiards. It clearly shows that the 
NNSD for the inverse x 2 distribution pi nvx 2 (y, s) agrees 
in a similar or, especially for small spacings, even better 
quality with the experimental ones as the others includ- 
ing the Brody (not shown) and the BR distribution. 
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TABLE II: Best-fit parameters for the experimental NNSD and resonance-strength distributions in the mixed billiards. The 
corresponding x 2 - vam es are given in brackets. 





Mushroom 








Limagon 






small 




large 




A = 0.125 




A = 0.150 


NNSD 
















x 2 


v = 2.95 (0.0189) 


IS 


— U.UO ^U.UUOl J 




— o en m w\ 


is 


— 9 (V\ (ft ftftOt7\ 


Inverse- x 2 


v = 0.00 (0.0042) 


V 


— i.Ol ^U.UUIOJ 


V 


— U.UU l^U.UUOl J 


is 


— ft nn ift nn9i "i 

— U.UU ^U.UUZIJ 


Log-normal 


v = 1.41 (0.0129) 


V 


n ctfi fa nnQ*^ 

— u.yo ^u.uuooj 


V 


1 90 /n m n/i ~\ 


V 


— 1.1/ ^u.uu to) 


Berry-Robnik 


q = 0.24 (0.0007) 


q 


= 0.87 (0.0027) 


q 


= 0.58 (0.0020) 


q 


= 0.62 (0.0034) 


Resonance strength distribution 














x 2 








IS 


= 2.66 (0.00142) 


IS 


= 3.47 (0.00035) 


Inverse- \ 2 








V 


= 0.50 (0.00030) 


V 


= 0.98 (0.00045) 


Gamma distribution 








IS 


= 0.76 (0.00120) 


IS 


= 0.81 (0.00131) 



TABLE III: Best-fit parameters obtained by extremizing the relative entropy for the experimental NNSDs the mixed billiards. 
The corresponding values of the relative entropy are given in brackets. 







Mushroom 






Limagon 








small 


large 




A = 0.125 


A = 0.150 


NNSD 














x 2 


V 


= 2.67 (-0.0920) v 


= 7.93 (-0.0232) 


V 


= 2.62 (-0.0909) v 


= 2.69 (-0.0827) 


Inverse- \ 2 


V 


= 0.01 (-0.0353) v 


= 3.14 (-0.0175) 


V 


= 0.01 (-0.0341) v 


= 0.01 (-0.0362) 


Log-normal 


V 


= 1.31 (-0.0568) v 


= 0.60 (-0.0206) 


V 


= 1.33 (-0.0533) v 


= 1.26 (-0.0537) 


Berry-Robnik 


q 


= 0.55 (-0.0289) q 


= 0.87 (-0.0203) 


q 


= 0.57 (-0.0293) q 


= 0.58 (-0.0378) 



The distribution p\ nvx 2 (is, s) (see Eq. (|16|) ) coincides 
with the Wigner distribution in the limit of is — > oo. As is 
decreases, the distribution evolves towards a well-defined 
limit, but this limiting case does not resemble the Pois- 
son distribution as one would expect. To demonstrate 
this behavior and give a feeling for the size of the tuning 
parameter is, we evaluate its value corresponding to the 
minimal deviation from a BR distribution with a given 
degree of chaoticity q. For this we define a measure 

f°° 2 

rflnv X 2 ,BR(^, q) = mill / [pi nvx 2(is,s) -pBB.(q, s)] ds. 

Jo 

(20) 

The distance di nvx 2 BR between the two distributions 
equals zero for q = 1 and v — * oo, where both distri- 
butions coincide with the Wigner distribution. Its value 
increases on departure from these parameter values, has a 
maximum value of 0.0030 for is = 1 and then decreases to 
a value of 0.0025 at is = 0. The measure <ij nvx 2 BR (z/, q) 
yields a relation between the BR and the superstatisti- 
cal parameter is, which is shown in Fig. [6l The figure 
suggests that j>i n vx 2 (^J s ) can om y describe the initial 
stage of the transition from chaos to regularity. One can 
find a BR distribution that agrees well with pi nvx 2 un- 
til an intermediate situation in which the NNSD corre- 
sponds to a Berry-Robnik parameter q = 0.83. The fail- 
ure of spacing distributions interpolating between those 
describing chaotic and regular systems in the limit of 



1.0 




0.6 I ■ 1 . 1 ■ 1 ■ 1 ■ 1 

2 4 6 8 10 

v 

FIG. 6: The relation between the tuning parameters q and v 
of the Berry-Robnik and superstatistical NNSDs, respectively. 



near integrability is common in different RMT descrip- 
tions based on generalized statistical mechanics, e.g. in 
[T^ . [20l [2ll. [2^. [23l [2^ ] . As mentioned above, the supersta- 
tistical random-matrix ensemble is base invariant. The 
Hamiltonian of an integrable system, on the other hand, 
by definition, has a well defined complete set of eigen- 
states, which constitutes a preferred basis. The RMT 
approach to mixed systems cannot depart far from the 
state of chaos without breaking base invariance. 
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Superstatistics 
Semi-Poisson 




FIG. 7: The superstatistical NNSD Pi nvx 2(0, s) compared 
with the semi-Poisson distribution Psp. 



The variance a 2 of the NNSD is often regarded as a 
one-parameter interpolation between chaos and order be- 
cause it monotonically increases from (4/V — 1) = 0.273 
for the Wigner distribution to 1 for the Poissonian. At 
v = 0, the superstatistical distribution in Eq. ([9]) becomes 



7T / 7T 

Plnv X 2 (°> s ) = ~ sK {2 s 



(21) 



and has a variance a 2 = (l6/ir 2 — l) = 0.621. This value 
is slightly larger than 0.5 which is exactly the variance of 
the semi-Poisson distribution 



psp(s) = 4sexp(-2s). 



(22) 



The semi-Poisson distribution was suggested to describe 
a narrow intermediate region between insulating and con- 
ducting regimes exemplified by the Anderson localization 
model [7l| . with the two limiting cases being described 
by Poisson and Wigner statistics, respectively. It was in- 
troduced to mimic new seemingly universal properties in 
certain classes of systems, in particular, being character- 
istics of the "critical quantum chaos" . Figure [7] compares 
the limiting superstatistical distribution, Pi nvx ^ (0, s), 
with the semi-Poisson distribution. We see from the fig- 
ure that the two distributions are quite similar in shape; 
the superstatistical distribution seems to just have passed 
the semi-Poissonian before reaching its final shape for 
v = 0. Therefore it is worth in the present context to 
use the semi-Poisson statistics as a reference distribu- 
tion marking the limit of validity of the base-invariant 
random-matrix description of mixed systems. 



V. RESONANCE STRENGTH 

In the preceding section we have seen that superstatis- 
tics yields a theoretical description of statistical proper- 
ties of the eigenvalues of a quantum system with mixed 
classical dynamics, which is of similar quality as existing 
ones, like e.g. the Brody and the BR ansatz. The great 



advantage of the method of superstatistics is, that we 
can similarly apply it to calculate the distribution of the 
eigenvector components of a Hamiltonian describing a 
mixed system in which the spectrum is composed of sub- 
spectra associated with levels following Poissonian and 
Wigner statistics (72[. For a chaotic system the squared 
eigenvector components follow the Porter-Thomas distri- 
bution 



PMt) = J—^\ 
V n y 



(23) 



where nis & constant parameter related to the mean 
value [6{|. In [72j superstatistics of the transition 
matrix elements was introduced by representing the 
transition-intensity distribution as a superposition of 
Porter-Thomas distributions with different values for the 
parameter 77. Similarly, in a superstatistical description 
of the squared eigenvector components of the Hamilton 
operator or equivalently the partial widths of the res- 
onances of a microwave resonator, the parameter 7/ in 
Eq. (j2"3"|) is no longer considered to be a constant but al- 
lowed to fluctuate according to a distribution f(rj). The 
superstatistical distribution of the squared eigenvector 
components is then given by 



f(v)\l^e-^drj. 



(24) 



The parameter distribution f(rj) was determined in [72j | 
using maximum-entropy arguments (33|. There, the re- 
sulting expression for the transition intensity distribution 
was fit to the distributions of the experimental reduced 
transition probabilities in 26 Al and 30 P nuclei [Z^, [z3] ■ It 
fits the data much better than a \ 2 distribution with v 
degrees of freedom, which is given in Eq. dHD with Vo = 1 
and has been proposed by Alhassid and Novoselsky (7f| 
and frequently used for describing the deviation of partial 
width distributions from the Porter-Thomas distribution. 
The superstatistical distributions obtained from Eq. ([24]) 



for a x j inverse x an d log-normal distribution f(rj) (see 
Eqs. (fg|- (fP0")0 . respectively, read 



-(u+l)/2 



(25) 



P 



(i+,)/2 (\/t{v + 2)) , 



and 



Plogniv^) 



\f2titv Jq y/fj 



: exp 



In I 2?7e 



-v 2 /4 



2v 2 



(27)" 

Unfortunately the latter integral could not be evaluated 
analytically. 



11 



In [76|, a new statistics, the resonance strength dis- 
tribution was introduced. For each resonance, a mea- 
surement of the transmission of microwave power from 
one antenna to another provides the product of the two 
partial widths related to the antenna "channels". This 
product has been named the strength of the resonance 
with respect to the transmission between the antenna 
channels. In Ref. [76| the strength distribution was in- 
vestigated experimentally for all three Limagon billiards. 
Altogether four antennas were attached to each of the mi- 
crowave billiards and the transmission spectra between 
pairs of antennas (a,b) were measured for all six pos- 
sible antenna combinations. An analytic expression for 
the strength distribution of a pair of partial widths with 
distributions P a and Ph is obtained as 



P(y)= / P a (QPb(tb)S(y - t a t h )dt a dt h . (28) 
Jo 

The distribution of the partial widths of a chaotic mi- 
crowave billiard is well described by the Porter-Thomas 
distribution Eq. (|23[) . and the corresponding resonance 
strength distribution is given as 



^GOE(y) = 



(29) 



More generally, if the partial width distribution is a 
X 2 distribution as proposed by Alhassid and Novoselsky, 
then Eqs. ([28]) and © aith rj = 1 yield 



Pa 



lh-Nc 



2 1 ~ v { Vy /y) v K Q { Vy /y) 

y r 2 (1//2) 



(30) 



The superstatistical resonance-strength distributions are 
obtained by substituting the corresponding partial- width 
distributions (Eqs. lj2"5 ]) -([2"7] ) ) into Eq. (fJH]). If the pa- 
rameter r\ has a log-normal distribution this leads to a 
double-integral which is not easy to calculate. For the 
case of a \ 2 distribution (see Eq. ([5])), the resonance- 
strength distribution is given by 



1 



(i/-2)7rV«F»(i//2) U"-2) 



II 



>+l)/2 

(31) 



r<22 

^22 



{y-2f 



1,1 



u+l 
2 • 



v+l 



where G™« 



ax, . 
h,. 



;bq 



is Meijer's G-function [67|, 



1771 . |78| |. On the other hand, if the parameter r\ has an 
inverse \ 2 distribution, one obtains 



Plnvx 2 ^^) 



47w 2 r 2 (1//2) 
1 



(32) 



/--f40 



16 



1 v v 
2' 2' 2 



Figure [5] compares the experimental resonance strength 
distributions for the two mixed Limagon billiards, ex- 




FIG. 8: The resonance strength distributions in the mixed 
Limagon billiards compared with the predictions of different 
superstatistics. The solid lines represent the x 2 , the dashed 
lines the inverse x 2 distribution, while the short-dashed line 
corresponds to the case in which the strengths have a \ 2 dis- 
tribution. 



pressed as functions of z = logio Hi with the correspond- 
ing distribution given in Eqs. ([30]), ((32]) and ((33j) . The 
best-fit values of the parameters are given in TabldlTl 
Both cases demonstrate the superiority of the supersta- 
tistical inverse \ 2 distribution. 



VI. SUMMARY 

Superstatistics has been applied to the study of a wide 
range of phenomena from turbulence to topics in econo- 
physic. RMT of the Gaussian ensembles is among these. 
In its application, the variance of the distribution of the 
matrix elements is chosen as a parameter whose distribu- 
tion is obtained by assuming a suitable form or applying 
the principle of maximum entropy. In this paper we use 
the time-series method to show that the spectra of mixed 
systems have two correlation scales as required for the 
validity of the superstatistical approach. The time-series 
analysis also shows that the best choice of the supersta- 
tistical parameter distribution for a mixed system is an 
inverse % 2 distribution. We computed the correspond- 
ing NNSD of the energy levels and compared it with 
the spectrum of two microwave resonators of mushroom- 
shaped boundaries and two of the family of Limagon 
billiards, which exhibit mixed regular-chaotic dynamics. 
The agreement is found to be similar to that with all the 
other considered well-established distributions including 
the celebrated BR distribution. The method of super- 
statics also provides a description of statistical proper- 
ties of the eigenfunctions of a system with mixed classi- 
cal dynamics. Thus, the resonance-strength distributions 
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for the Limacon billiards could be analyzed. The agrree- 
ment with the experimental resonance-strength distribu- 
tions is better than with that derived from the Alhassid- 
Novoselsky x 2 distribution of transition intensities. 
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